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1  INTRODUCTION 


1.1  Program  Overview 

The  purpose  of  the  FY92/93  effort1  was  to  investigate  and,  where  appropriate,  develop  and  de¬ 
liver  software  tools  to  improve  the  productivity  and  quality  of  oceanographic  data  interpretation. 
Specifically,  this  effort  comprised  three  tasks: 

•  Task  1 :  Research  and  Develop  Methods  tor  Improving  the  Performance  of  Automated  Outlier 
Detection.  Based  on  algorithms  and  software  developed  by  SRI  under  a  previous  DMA  pro¬ 
gram,2  several  enhancements  were  carried  out. 

•  Task  2:  Implement,  Document,  and  Deliver  the  Improved  Outlier  Detection  Algorithms.  Addi¬ 
tional  modifications  were  made  to  support  generic  multibeam  format  (GMF)  data  sets.  The 
upgraded  software  was  delivered  to  the  hydrographic  and  bathymetric  groups  at  the  Naval 
Oceanographic  Office  (NAVOCEANO)  and  NOAA,  and  will  be  delivered  to  the  Defense 
Mapping  Agency  (DMA)  concurrent  with  this  final  report.  These  improvements  are  summa¬ 
rized  in  Section  2. 

•  Task  3:  Identify  and  Evaluate  More  Advanced  Oceanographic  Data  Analysis  Tools.  A  study 
was  made  to  determine  the  data  automation  needs  ot  the  oceanographic  community,  and 
then  compare  these  needs  to  data  analysis  tools  already  developed  for  other  scientific 
application  areas.  The  survey  of  data  automation  needs  was  carried  out  with  the  assistance 
of  NAVOCEANO,  DMA,  and  NOAA.  Our  research  indicated  that  a  significant  body  of  rele¬ 
vant  techniques  has  indeed  been  developed  for  applications  such  as  robotics,  image  under¬ 
standing,  and  various  statistical  analysis  applications.  A  summary  of  our  results  was  pre¬ 
sented  in  our  February  1993  technical  research  memorandum,  “Oceanographic  Data 
Processing  Feasibility  Study”  1 1]. 

1.2  Background 

In  a  previous  (FY90/91)  effort  sponsored  under  the  DMA  Defense  Hydrographic  Initiative 
(DHI),  SRI  was  tasked  to  develop  of  a  general-purpose  outlier  detection  algorithm  for  use  with 
multibeam  sonar  systems.  Depth  outliers  are  a  particularly  common  problem  in  deep-water  multi¬ 
beam  sonar  systems,  and  must  be  eliminated  from  the  data  prior  to  seafloor  map  production.  The 
time  required  to  manually  locate  and  remove  outliers  represents  a  significant  portion  of  the  overall 
data  production  time  and  labor.  Furthermore,  missed  outliers  can  lead  to  substantial  errors  in  the 
final  map  product. 

The  result  of  this  earlier  effort  was  the  development  and  deployment  of  a  software  tool  that 
automates  outlier  recognition  and  deletion.  Two  diverse  approaches  to  outlier  detection  were 
developed:  the  Graduated  Non-Convexity  (GNC)  algorithm,  adapted  from  computer  vision 


1.  1992/93  contract  managed  by  the  Office  of  Naval  Research,  Contract  N00014-92-C-0015. 

2.  1990/91  contract  also  managed  by  the  Office  of  Naval  Research,  Contract  N00014-90-C-0132. 
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applications;  and  the  Robust  Linear  Prediction  (RLP)  algorithm,  adapted  from  ihe  field  of  digital 
image  restoration.  Both  algorithms  were  consolidated  in  a  single  interactive  software  package  that 
included  an  interface  to  several  bathymetry  data  formats,  parameter  controls,  and  graphic  dis¬ 
plays.  Prototype  software  releases  were  delivered  to  the  Naval  Oceanographic  Office 
(NAVOCEANO),  NOAA,  and  DMA. 

The  performance  of  this  system  was  investigated  in  experiments  on  U.S.  Navy  SASS,  Sea- 
beam,  and  Krupp  Atlas  Hydrosweep  bathymetry  data  sets.  An  extensive  analysis  was  carried  out 
on  a  collection  of  representative  SASS  data  sets  by  a  NAVOCEANO/GBA  team.  In  the  NAV¬ 
OCEANO  experiment,  the  outliers  located  by  the  algorithms  were  compared  with  those  identified 
by  experienced  Navy  bathymetrists.  The  automatic  outlier  detection  system  located  66%  -  85%  of 
valid  outliers  (where  the  performance  varied  somewhat  depending  on  the  nature  of  the  test  data 
sets),  with  false-detection  rates  of  less  than  0.03%  [21.  These  results  were  considered  encouraging 
by  the  oceanographic  community.  The  low  false-detection  rate  is  especially  important,  because 
this  ioss  rate  is  small  in  comparison  to  typical  sonar  data  dropout  rates  and  is  believed  not  to  ma¬ 
terially  alter  the  final  bathymetry  product.  As  a  consequence  of  these  results,  NAVOCEANO  has 
integrated  the  outlier  detection  algorithms  into  their  bathymetric  data  production  process. 

These  tests  also  identified  high-priority  performance  improvements  for  the  automated  outlier 
detection  algorithms.  Specific  suggestions  were  made  by  groups  in  the  oceanographic  community 
(as  reported  in  Arnold  and  Shaw  [2]),  and  these  form  the  basis  for  the  algorithm  research  and  de¬ 
velopment  carried  out  under  Task  1  of  the  current  program.  An  additional  requirement  to  support 
NAVOCEANO  hydrographic  data  production  was  recognized,  and  set  forth  as  Task  2.  Finally,  a 
requirement  to  automate  other  data  interpretation  tasks,  as  called  for  by  the  DHI  program,  is  the 
basis  for  the  feasibility  study  carried  out  under  Task  3. 

1.3  Acknowledgments 

As  in  our  previous  effort,  SRI  received  technical  assistance  that  proved  to  be  absolutely  essen¬ 
tial  to  our  work.  We  would  like  to  thank  the  staff  of  NAVOCEANO,  and  in  particular  Rich  Sandy, 
who  took  the  time  to  tutor  us  in  the  special  needs  of  bathymetric  and  hydrographic  data  produc¬ 
tion;  and  Dr.  Tom  Davis,  now  of  S  AIC,  for  his  highly  relevant  comments  and  continued  support  of 
automation  effort'. 
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2  TASK  SUMMARIES 

2.1  Task  1 :  Outlier  Detection  Algorithm  Improvements 

Two  fundamentally  different  outlier  detection  algorithms  were  considered  under  this  program. 
The  first  technique,  the  graduated  non-convexity  (GNC)  algorithm,  reconstructs  the  ocean  floor 
surface  based  on  a  combination  of  the  sonar  depth  estimates  and  constraints  lor  surface  smooth¬ 
ness  and  discontinuity.  Sonar  sample  points  that  significantly  differ  from  the  re-estimated  surface 
are  identified  as  potential  outliers.  The  second  approach,  the  robust  linear  prediction  (RLP)  algo¬ 
rithm,  adaptively  develops  statistical  model  parameters  from  the  data  set.  The  expected  depth  of 
each  sample  point  is  predicted  from  the  values  and  statistics  of  the  neighboring  depth  samples.  If 
the  actual  depth  varies  significantly  from  the  predicted  value,  the  sample  is  identified  as  an  outlier 
and  the  region  is  reviewed  with  the  anomalous  sample  removed. 

2.1.1  GNC  Algorithm  improvements 

The  outlier  detection  properties  of  the  GNC  are  a  side  effect  of  the  algorithm’s  ability  to  de¬ 
velop  a  surface  with  discontinuities:  the  system  has  been  adapted  to  recognize  very  small  and 
hence  physically  impossible  discontinuities  as  outliers.  The  computer  vision  model,  unfortu¬ 
nately,  does  not  accurately  reflect  the  true  characteristics  of  the  ocean  floor,  which  is  typically  a 
continuous  (albeit  sometimes  steep)  surface  that  rarely  contains  occlusions. 

A  second  point  is  that  the  original  (FY90/91 )  GNC  algorithm  was  based  on  a  first-order  statis¬ 
tical  model  for  the  seafloor.  In  a  first-order  model,  the  depth  differences  between  adjacent  points 
on  the  seafloor  are  modeled  as  a  zero-mean,  Gaussian  random  process.  Review  of  statistical 
studies  by  Fox  [3]  indicated  that  this  model  is  accurate  in  relatively  flat  seafloor  regimes,  but  that 
a  second-order  model  is  more  appropriate  in  regions  with  strong  topographic  features  (Figure  1). 
This  finding  corroborates  our  1991  experimental  evidence  that  the  first-order  model  had  difficulty 
in  steeply  sloped  regions  [2]. 

From  the  above,  it  was  determined  that  the  surface  energy  function  used  by  the  GNC  should 
be  reformulated  to 

1.  Emphasize  singularities  (outliers). 

2.  Incorporate  a  second-order  energy  model. 

Background  on  GNC  Algorithm.  The  GNC  algorithm  is  a  visible-surface  reconstruction 
technique  originally  developed  for  computer  vision  and  three-dimensional  mapping  systems.  In  a 
typical  GNC  application,  such  as  the  simple  example  of  Figure  2,  an  active  range-finder  sensor 
(sonar  or  an  IR  ranging  device)  develops  a  collection  of  range  estimates  over  some  viewable 
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FIGURE  1  SPECTRAL  MODEL  COMPARISON 


Field  of  view 


(a)  Example  Computer  Vision  Application 


(b)  Corresponding  Raw  Depth  Map 


(c)  Estimated  Visible  Surface 

FIGURE  2  VISIBLE-SURFACE  RECONSTRUCTION 
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region  of  interest.  After  the  sensor  performs  a  scan,  a  potentially  noisy  discrete  depth  map  is 
created.  The  resulting  visible  surface  consists  of  piecewise  continuous  regions  formed  by  object 
surfaces,  with  discontinuities  (i.e.,  significant  jumps  in  depth)  at  the  object  boundaries. 

The  GNC  algorithm  is  used  to  estimate  the  object  edge  locations,  as  w'ell  as  to  reduce  the 
noise  and  thereby  improve  the  shape  and  position  estimates  of  object  surfaces.  A  combination  of 
factors  is  used  in  the  surface  reconstruction  process,  including  constraints  for  smoothness  (e.g., 
prior  knowledge  on  the  shapes  of  potential  objects),  depth  discontinuities,  and  reliabi'ity  of  the 
measurement  data.  These  competing  criteria  are  quantified  in  terms  of  a  surface  energy  function, 
where  the  minimum-energy  surface  represents  the  most  probable  solution.  The  primary  technical 
contribution  of  the  GNC  algorithm  is  the  nonlinear  optimization  procedure  by  which  it  locates  the 
minimum-energy  configuration. 

Surface  Energy  Models.  GNC  constuins  the  estimated  surface  by  applying  simultaneous  con¬ 
straints  for  surface  smoothness,  continuity,  and  adherence  to  the  raw  sonar  data,  as  quantified 
through  the  function  E  ( u ,  d)  : 

E(u,d)  -  D  (u,d)  +  S  (u)  +  P  (u)  (1) 

where 

D  («,  d)  measures  adherence  of  the  smoothed  surface  estimates  to  the  data 
S  (u)  is  a  surface  smoothness  function 
P  (u)  is  a  surface  and/or  derivative  discontinuity  function. 

The  energy  model  used  in  the  FY90/91  GNC  implementation  is  based  on  a  first-order,  or 
“weak  membrane”  model.  The  weak-membrane  analogy  is  a  thin  rubber  sheet  that  attempts  to 
minimize  the  vertical  slope,  but  tears  if  stretched  beyond  some  strength  limit.  The  weak-mem¬ 
brane  energy  function  is  calculated  as  a  sum  of  energy  contributions  over  the  entire  surface  under¬ 
going  reconstruction: 


where 


Ej  (u,  d) 


X  { (w9y  “  dcj)  2  +  ( 1  ~  Li,j)  x  +  “P  +  aLi,j 1 

i,j  e  S 


(2) 


£j  ( u ,  d)  is  the  first-order  surface  energy  function 

U-  j  and  d-  ■  are  discrete  samples  from  the  reconstructed  surface  and  sonar  measurement  set, 
respectively,  at  each  point  { i,  j }  on  the  reconstructed  surface,  S 

\  is  a  first-order  smoothness  parameter. 


o 


a  is  the  energy  associated  with  surface  discontinuities  (analogous  to  the  membrane 
“strength”). 

ux  and  uy  are  the  x  and  y  spatial  first  derivatives  at  each  point  { i,  j }  . 

L-  ■  is  a  discontinuity  indicator,  such  that: 

_  1  at  discontinuities 

[  0  elsewhere 

The  variable  Li  ■  enables  surface  discontinuities  that  are  key  to  this  application;  at  points 
where  Li  .  -  1 ,  the  smoothness  constraint  is  ignored  and  a  fixed  energy,  a,  is  assigned.  At  all 
other  points,  the  smoothness  energy',  A2  (u2  +  ul) ,  is  assessed.  The  relative  value  of  the  smooth- 
ness  parameter,  X ,  and  the  discontinuity  energy,  a,  determine  the  algorithm’s  proclivity  toward 
discontinuities. 

A  second-order,  or  tin  plate,  model  was  developed  under  the  1992/93  eifort.  The  second- 
order  approach  provides  surface-fitting  behavior  like  that  of  a  stiff  but  deformable  plate.  As  in¬ 
ferred  from  the  thin-plate  analogy,  this  model  favors  surfaces  with  constant  slopes,  while  resisting 
strong  surface  curvature.  The  surface  is  allowed  to  break,  however,  causing  discontinuities  in  the 
surface  and  its  first  derivative.  The  thin-plate-model  energy  function  used  in  this  effort  is  given  as 

EZ  (u’  =  £  1  {uiJ  ~  di'P  "  +  ^  ( 1  ~  LiJ  (i4  +  ujy>  +  aLij\ 

ijeS  (3) 

where 


E2  ( u ,  d)  is  the  second-order  surface  energy  function 

p.  is  the  secord  order  smoothness  parameter 
a  is  the  energy  associated  with  surface  discontinuities 
u  and  m  are  the  second-order  spatial  derivatives  at  points  { i,  j} 
Li  j  is  the  discontinuity  indicator. 


Surface  Energy  Minimization.  Reliable  minimization  of  the  energy  function,  E  («,  d) ,  requires 
a  complete  search  over  all  possible  reconstructed  ocean  depths,  «•  ,  for  all  points  on  the  surface. 
In  most  applications,  and  in  particular  bathymetric,  such  searches  are  intractable  because  of  the 
extremely  large  number  of  possible  surfaces.  Alternatively,  if  the  locations  of  the  discontinuities. 
L-  ,  are  fixed,  the  minimization  can  be  posed  as  a  common  linear  least-squares  problem.  In  this 


3.  Other  second-order  energy  functions  arc  also  possible;  for  example,  a  reasonable  function  might  account  for  both 
first-  and  second-order  smoothness,  and/or  rely  on  the  derivative  crossierms  uiy  and  u>i. 
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case,  a  least-squares  solution  to  u-  /  would  be  found  for  all  possible  collections  discontinuities, 
and  the  discontinuity  set  providing  the  minimum-energy  surface  would  be  chosen.  This  may  be  a 
plausible  approach  if  the  number  of  expected  ocean  floor  discontinuities  is  small.  The  processing 
requirements  quickly  grow'  unmanageable,  however,  because  the  number  of  feasible  outlier  con¬ 
figurations  grows  combinatorially  with  the  number  of  possible  outliers. 

Blake  and  Zisserman  proposed  the  GNC  as  a  practical  solution  for  problems  of  this  kind  |4], 
The  primary  contribution  of  the  GNC  algorithm  is  the  unique  method  by  which  it  carries  out  the 
energy  minimization  process.  Rather  than  minimize  the  energy  function  directly  (e.g.,  through  an 
exhaustive  search),  the  algorithm  forms  a  bounded  approximating  function  with  special  convexity 
properties.  Once  the  minimum  of  the  approximating  (convex)  function  is  located,  a  more  accurate 
approximating  function  is  created  and  the  minimization  process  is  repeated.  The  accuracy  of  the 
energy-minimum  estimate  is  refined  as  the  approximating  function  approaches  the  true  energy 
function. 

For  notational  convenience,  the  smoothness  and  discontinuity  contributions  of  the  energy 
function  are  usually  consolidated,  forming  the  energy  function 

F  =  D  +  L[g(uxx)  +g(u>7)],  (4) 

where  D  is  the  measurement  adherence  function,  and  g  (r)  contains  both  the  smoothness  and  dis¬ 
continuity  functions.  In  the  case  of  the  second-order  model,  g  (i)  is  given  by 

8(0  =P2( 1  -Li  j)t2  +  aLi  j .  (5) 

Blake  and  Zisserman  have  shown  { 4 ]  that  the  best-fit  convex  approximation  to  g  ( i )  is 
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for  (  r  >  r) 
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r  =  a  (  -  +  — 
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q  = 
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P  *■  4'  m 

and  t  is  the  second  derivative  about  the  depth  sample.  The  function  g  ^  takes  any  of  three  forms, 

depending  on  the  value  of  the  argument.  For  small  depth  differences,  the  magnitude  of  /  is  small 
and  the  smoothness  energy  criteria,  p4;2,  applies.  Large  depth  differences  represent  discontinui¬ 
ties,  and  the  discontinuity  energy,  a,  is  assessed.  In  both  the  small-  and  large-difference  cases,  the 
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approximating  function  is  exact.  At  intermediate  depth  differences,  a  quadratic  (and  hence 
convex)  approximation  is  used. 

The  parameter  p  controls  the  energy  function  approximation.  The  algorithm  locates  the  solu¬ 
tion  iteratively,  beginning  with  p  =  1 ,  such  that  g  1 1  ’  =  g*  and  the  approximating  energy  func¬ 
tion  is  globally  convex.  After  a  solution  for  g  {  ^  is  located,  p  is  decreased  and  a  new  (more  accu¬ 
rate)  solution  is  found.  As  the  algorithm  proceeds,  p  is  slowly  decreased  toward  zero  and  a  refined 
estimate  is  calculated  at  each  point.  Note  that  at  p  =  0,  the  approximating  energy  function  is  pre¬ 
cise,  and  the  true  eneigy  minimum  is  located.  The  algorithm’s  operation  is  summarized  in  Table  1. 

_ Table  1.  SECOND-ORDER  GNC  ALGORITHM  SUMMARY _ 

Initialize: 

for  i  =  1, ...  N; ,  j  =  1, ...  Nj 


Reconstruct: 

for  p  =  1,1/2, 1/4, ... 

for  n  =  1, ...  N  (iterauons  on  the  reconstructed  image) 
for  i  =  1, ...  Nj  j  *  1  ...  N3 


<j 1 


(p) 


where  u**  and  u\y-  are  the  second  derivatives  with  respect  to  depth  in  the  x  and  y  dimensions,  re¬ 
spectively,  of  the  reconstructed  seafloor  surface,  u,  at  iteration  interval  n ;  dt }  is  the  raw  depth 
sample  at  position  {i,j} ;  and  y  is  the  adaptation  gain  constant. 

Empirical  Comparison  of  First-  and  Second-Order  Models.  Seafloor  reconstruction  results  ob¬ 
tained  with  the  first-  and  second-order  GNC  algorithms  are  given  in  Figure  3.  The  test  data  set 
was  excerpted  from  a  1990  University  of  Rhode  Island  Hydrosweep  cruise  of  the  North  Atlantic 
near  64°latitude,  5°  longitude.  As  shown  in  the  three-dimensional  plot  of  Figure  3a.  the  selected 
data  set  contains  a  significant  topographic  feature.  Results  from  the  first-order  GNC  algorithm 
(Eq.  [3])  are  shown  in  Figure  3b.  Note  that  the  first-order  reconstructed  surface  has  linear  tears 
following  across  the  steeply  sloped  regions,  causing  the  reconstructed  surface  to  differ  greatly 
from  the  underlying  data.  True  outliers  in  the  steep  regions  are  masked  by  these  tears,  thereby  re¬ 
ducing  the  outlier  detection  effectiveness. 
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The  second-order,  or  thin-plate  reconstruction  (Eq.  (4|)  results  are  shown  in  Figure  3c.  Be¬ 
cause  the  second-order  model  is  a  better  match  to  this  data  set,  the  reconstructed  surface  is  contin¬ 
uous  through  steep  regions,  yet  identifies  breaks  in  the  surface  to  accommodate  outliers.  Several 
outliers  are  present,  however,  in  the  form  of  occasional  spikes  in  the  data. 

2.1.2  RLP  Algorithm  improvements 

The  RLP  algorithm  operates  by 

1.  Estimating  autoregressive  parameters  for  some  window  of  data. 

2.  Predicting  depth  at  a  beam-ping  bin  based  on  the  estimated  parameters  and  some  neighbor¬ 
hood  of  surrounding  bins. 

3.  Comparing  the  predicted  depth  to  the  measured  depth  at  that  bin  and  making  some  decision. 

Rather  than  being  arbitrary,  the  shape  of  the  prediction  neighborhood  is  governed  by  the  math¬ 
ematics  of  the  parameter  estimation  procedure.  To  simplify  this  procedure,  we  followed  standard, 
image-processing  procedure  and  employed  a  nonsymmetric  prediction  neighborhood.  Unfortu¬ 
nately,  this  choice  has  unwanted  side  effects,  including: 

•  Probability  of  outlier  detection  depends  on  the  direction  of  filtering,  whether  from  early 
times  to  later  times,  or  vice  versa. 

•  Because  of  the  one-sided  prediction,  the  RLP  tends  to  produce  false  alarms  at  transition 
zones. 

Improvements  to  the  RLP  algorithm  were  designed  to  overcome  the  directionality  of  the  non- 
symmetric-half-plane  support  region  and  to  lower  the  false-alarm  rate  on  steeply  dipping  regions. 

Background  on  RLP  Algorithm.  The  RLP  algorithm  characterizes  outliers  as  samples  of  a  non- 
Gaussian  random  process.  The  tools  of  robust  estimation  theory  allow  us  to  detect  these  samples 
without  assuming  a  parametric  form  of  the  outlier  probability  density.  The  RLP  algorithm  oper¬ 
ates  adaptively  by  robustly  estimating  the  parameters  of  an  autoregressive  (AR)  ocean  bottom 
process  within  some  small  window,  and  then  using  these  parameters  to  predict  depths  within  the 
same  window.  The  estimation  window  steps  across  a  few  beams  and  pings  at  a  time,  as  shown  in 
Figure  4.  If  a  measured  depth  exceeds  some  threshold,  it  is  replaced  by  its  predicted  value  and  the 
estimation  process  repeats.  Outliers  are  declared  when  the  difference  between  the  restored  surface 
and  the  original  surface  exceeds  some  threshold.  This  approach  offers  several  benefits: 

•  The  ocean  bottom  is  modeled  accurately,  even  in  the  presence  of  large  sensor  errors. 

•  What  is  known  (the  ocean  bottorr  statistics)  is  exploited,  and  what  is  unknown  (the  outlier 
statistics)  is  removed. 
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FIGURE  4  MOVEMENT  OF  ESTIMATION  WINDOW  ACROSS  DATA  SET  IN  OVER¬ 
LAPPING  STEPS 

•  The  process  adapts  to  the  local  ocean  bottom  statistics  without  reliance  on  a  global  surface 
estimate. 

•  Ocean  bottom  discontinuities  are  modeled  and  not  reported  as  outliers. 

•  Allowance  is  made  for  random  sampling  errors  (beam-ping  format  data). 

Robust  Bathymetric  Error  Detection.  The  RLP  algorithm  was  adapted  from  an  image  restora¬ 
tion  algorithm  reported  by  Kashyap  and  Eom  [5].  In  this  technique,  a  robust  M-estimator  is  used 
to  estimate  the  parameters  of  an  autoregressive  image  model.  These  parameters  are  then  used  in  a 
“data  cleaning”  step,  where  outliers  are  removed.  The  AR  parameters  are  again  estimated  and  the 
data  cleaned,  and  so  on,  until  no  more  outliers  are  found. 

The  Kashyap  algorithm  assumes  a  first-order,  nonsymmetric  half-plane  (NSHP)  image  model. 
This  model  is  the  two-dimensional  equivalent  of  a  causal  model  in  time  series  analysis  (Figure  5). 
In  the  NSHP  model,  an  intensity,  y  (i,j) ,  is  computed  as  the  weighted  sum  of  the  three  “neigh¬ 
bors”  shown  in  Figure  5,  plus  some  random  innovation.  If  0^.  is  the  weight  given  to  Lhe  k‘h 
neighbor  of  point  (i,  j) ,  this  relation  may  be  written 

yVJ)  =  8jy(i  -  l,y)  +02y(i~  \,j-  1)  +03y  (ij-  1)  +  0(1,7) , 

where  x>  ( i,j )  is  a  zero-mean  Gaussian  random  Field  with  variance  o*.  In  vector  notation,  this 
would  be  written 


yU,j)  =  +\)(i,;). 


(7) 
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where  z  (i,j)  is  a  3x1  vector  of  NSHP  neighbor  values  surrounding  point  ( i,j)  .  The  image  is 
further  assumed  to  be  contaminated  by  some  unknown  outlier  process,  £  (/,/') ,  whose  probability 
distribution  is  unknown.  Thus,  the  observed  image,  y0(ij) ,  may  be  expressed  as 

y0(ij)  =  Qll(i,j)  + 1)  (/,/)  +C(U)  • 

Given  this  observation  model,  a  robust  M-estimator  for  the  parameters  0  and  the  variance  is 
found  by  minimizing  the  nonlinear  expression 


=  Sp 


The  data-cleaning  procedure  takes  the  process  further  by  replacing  outliers  in  each  step  with 
the  value  predicted  by  the  robustly  estimated  parameters  0.  Upon  convergence,  the  window  is 
shifted  over  a  few  rows  or  columns  and  the  process  repeated.  In  this  way,  the  entire  image  is  re¬ 
stored,  one  window  at  a  time. 

RLP  Modifications.  Throughout  development  of  the  RLP  algorithm,  we  have  relied  on  a  first- 
order  nonsymmetric  half-plane  image  model,  which  has  proven  sufficiently  flexible  to  model  a 
variety  of  ocean  bottom  processes.  It  also  has  the  advantage  of  being  recursive,  so  that  missing 
data  points  may  be  estimated  as  the  filter  progresses.  The  NSHP  model  also  has  several  disadvan¬ 
tages.  It  has  strong  directional  preferences,  and  only  detects  errors  on  the  leading  edge  of  a  group 
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•  •  o  oj«  • 

•  •  of®  •  • 


FIGURE  5  NSHP  MODEL:  UNFILLED  POINTS  COMPOSE  THE 
NEIGHBORHOOD  OF  THE  CIRCLED  POINT. 
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of  error  points.  This  directionality  problem  may  be  solved  through  the  use  of  noncausal  error  fil¬ 
ters. 

A  more  serious  drawback  of  the  NSHP  model  is  that  it  does  not  refleci  the  nature  of  the  multi¬ 
beam  bathymetry  process.  In  other  words,  the  elements  of  the  vector  z  ( i,j)  in  Eq.  (7)  do  not  cor¬ 
respond  to  beams,  pings,  or  any  physical  entity  other  than  an  arbitrary  neighborhood  of  points. 
This  lack  of  correspondence  makes  impossible  the  incorporation  of  a  priori  information,  such  as 
beam  or  ping  confidence  measures.  The  utility  of  the  RLP  has  been  improved  by  reformulating 
the  surface  model  to  be  oriented  specifically  to  multibeam  bathymetry.  In  this  algorithm,  which 
we  refer  to  as  the  vector  robust  linear  prediction  (VRLP)  algorithm,  beam  data  can  be  expressed 

as  a  vector,  Markov  random  process  by  replacing  Eq.  (2)  with 

N 

p{t)  =  £  Qpit-i)  +v(r) 

i  =  i 

where  p  (t)  is  a  length  n  vector  of  measurements  from  a  ping  at  time  t,  one  component  for  each 
beam.  The  parameters  of  this  model  are  contained  in  @.,  an  n  x  n  matrix,  and  v  (t)  is  a  vector 
random  process  with  covariance  matrix  C.  The  estimated  covariance  matrix  C  for  such  a  model 
now  conveys  much  more  information  than  the  covariance,  a,  of  the  scalar  NSHP  model  used  pre¬ 
viously.  For  example,  the  diagonal  of  C  contains  estimates  of  the  variance  of  individual  beams.  If 
a  single  element  along  the  diagonal  grows  much  larger  than  the  other  elements,  an  operator  can  be 
alerted  to  the  presence  of  a  bad  beam. 

Rather  than  directly  copy  the  Kashyap  algorithm,  we  have  adopted  a  vector  version  of  the  ro¬ 
bust  filter-cleaner  approach  of  Martin  and  Thomson  [6].  This  algorithm  is  identical  to  the 
Kashyap  approach  with  the  exception  that  simple  least-squares  estimates  (which  we  obtain  by 
solving  vector  version  of  the  Yule-Walker  equations  [7])  of  the  parameters  and  covariance  are 
computed  at  each  iteration  rather  than  taking  a  nonlinear  Gauss-Newton  step.  The  data  cleaning 
proceeds  as  before,  with  certain  modifications  for  the  vector  quantities  involved. 

To  improve  false-alarm  rates  on  steep  slopes,  we  have  further  modified  the  algorithm  to  per¬ 
form  smoothing  rather  than  filtering.  In  other  words,  two  models,  one  forward  and  one  backward, 
are  assumed,  and  a  weighted  sum  of  each  is  used.  This  is  shown  graphically  in  Figure  6.  The 
models  may  be  written 

Nt 

pf.  =  Y  &  +  ^ 

‘  j  £->  *,  o  -  0  j 

i  =  1 
and 
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pb  =  y  ©f  lPbr  ..  +df. 

W  Zj  i.  fc'-  0  +  0  7 

i  =  1 

where  t3^  and  are  vector  random  processes  with  covariance  matrices  C/  and  C*,  respectively. 
The  complete  VRLP  algorithm  is  given  in  Table  2. 


_ Table  2.  VECTOR  ROBUST  LINEAR  PREDICTION  ALGORITHM _ 

1.  k  =  0 

2.  Until  no  change  in  0  or  c  , 

do 

2. a  U^ing  data  within  the  estimation  window,  compute  0^  k ,  0*  ^  *,  and 

C,-  r  i  =  1 ,  ...,N  ,  by  solving  the  forward -backward,  vector,  Yule-Walker 
equations. 

2.b  Perform  data  cleaning: 

For  each  ping  ;n  the  estimation  window, 

do 


1  =  1 
N 


i  =  1 
N 


'  •  L  v  i  =  1 


i  =  1 

2.b.iii  Pjk+i  =  «7‘+CS1)-1(C7‘^ttl+CJ1^^1) 

end 


end 


2.c  it  =  *  +  1 


RLP  Results.  A  comparison  of  the  RLP  and  VRLP  algorithms  is  shown  in  Figure  7.  Figure  7a 
shows  a  raw  data  swath  with  a  variety  of  errors.  Figure  7b  shows  the  results  of  applying  the  RLP 
algorithm  to  these  data.  Note  that  several  false  alarms  are  found  on  the  steeply  sloping  region,  and 
that  some  outliers  are  missed  on  the  edges  and  in  the  vicinity  of  a  bad  beam.  Figure  7c  shows  the 
results  of  applying  the  VRLP  algorithm.  Note  that  the  number  of  false  alarms  is  reduced,  and  that 
more  outliers  are  detected  on  the  edges  and  in  the  middle. 


15 


Estimate  local 


FIGURE  6  BIDIRECTIONAL  VRLP  ALGORITHM 


2.2  Task  2:  Outlier  Detection  Software  Implementation  and  Delivery 

2.2.1  Background 

The  purpose  of  the  clean  program  is  to  automatically  analyze  raw  depth  files  produced  by 
bathymetric  and  hydrographic  sensors,  and  tag  the  isolated  outlier  depth  estimates  in  the  data.  The 
application  environment  for  clean  is  summarized  in  Figure  8,  where  clean  represents  the  first  step 
in  a  four-step  bathymetric  chart  production  process.  In  the  second  step  of  Figure  8,  the  positional 
information  attached  to  the  depth  estimates  is  adjusted  to  account  for  navigational  errors.  The 
third  step  calculates  and  assigns  statistical  information  to  the  data,  such  as  data  reliability  indica¬ 
tors  and  expected  errors  (in  position  and/or  depth).  The  fourth  and  final  step  combines  the  indi¬ 
vidual  swaths  of  depth  estimates  to  form  a  fully  populated,  two-dimensional  depth  grid. 

Clean  can  be  applied  either  as  a  stand-alone  prefilter  of  sensor  data,  or  used  in  conjunction 
with  the  NAVOCEANO  swathedt  contour  display  program  to  provide  a  semiautomated  oceano¬ 
graphic  data  quality  analysis  and  editing  environment.  The  auxiliary  program  control  provides  an 
X-Windows  interface  to  the  clean  processing  parameters. 

2.2.2  System  Operation 

A  top-level  system  control  and  data  flow  diagram,  drawn  from  the  operator’s  perspective,  is 
shown  in  Figure  9.  Input  to  the  system  consists  of  raw  sensor  depth  files  in  one  of  the  following 
formats: 
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FIGURE  8  APPLICATION  ENVIRONMENT:  BATHYMETRIC/ 
HYDROGRAPHIC  CHART  PRODUCTION 


•  Compressed  SASS  and  Seabeam 

•  Uncompressed  Seabeam 

•  Hydrosweep 

•  Generic  multibeam  format  (GMF). 

Clean  optionally  produces  two  types  of  output  files: 

•  Filtered  Time  Anomaly  File  (FTAF) 

•  Edited  sensor  depth  file. 

The  FTAF  is  an  ASCII  file  containing  the  times  and  beam  numbers  of  detected  outliers.  Other 
software,  including  the  swathedt  contouring  program,  makes  use  of  the  FTAF  to  graphically 
highlight  oudiers.  The  edited  sensor  data  file  has  identical  format  to  that  of  the  input  file,  except 
that  the  automatically  detected  outlier  depths  have  been  marked. 
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FIGURE  9  CLEAN  SOFTWARE  BLOCK  DIAGRAM 


Clean  is  a  batch  program;  once  started,  it  does  not  require  operator  interaction.  At  the  start  of 
execution,  clean  is  assigned  processing  parameters  that  include; 

•  Input  file  name  (required) 

•  FTAF  file  name  (optional;  if  omitted,  no  FTAF  created) 

•  Range  of  pings  and  beams  to  process  ‘  optional;  if  omitted,  all  pings  &  beams  processed) 

•  Outlier  detection  algorithm  parameters  (optional;  default  values  used  otherwise). 

Processing  continues  until  either  the  operator  aborts  the  run  (by  explicitly  killing  the  process 
via  an  operating  system  command,  e.g.,  the  Unix  kill  command),  or  until  the  specified  range  of 
pings  is  completely  processed. 

The  auxiliary  program  control ,  an  X-windows  interface,  enables  more  convenient  viewing  and 
editing  of  the  parameters  than  is  possible  through  standard  Unix  command-line  parameter 
passing.  In  addition,  control  allows  entire  sets  of  parameters  to  be  saved  and  recalled  using  oper¬ 
ator-selected  names.  The  parameter  values  selected  in  the  control  window  are  communicated  to 
clean  via  a  parameter  database  file  that  is  shared  between  the  two  programs.  An  example  control 
screen  image  is  shown  in  Figure  10. 

2.2.3  Computing  Requirements 

Clean  was  designed  to  operate  on  any  Unix  computing  system,  and  to  be  readily  portable  (that 
is,  transferable  with  minimum  effort)  to  other  operating  systems. Because  clean  is  a  computation¬ 
ally  intensive  program,  a  computer  with  fast  floating-point  performance  is  recommended.  Clean 
does  not  require  interactive  or  graphics  support  software.  The  auxiliary  program  control,  however, 
requires  a  graphics  display  terminal  with  keyboard  and  mouse,  along  with  an  X-Windows  server 
and  window  manager.  The  details  of  control's  appearance  depend  on  the  window  manager  and  the 
local  shareable  libraries  supporting  X.  Control  was  implemented  under  the  MOTIF  widget  set. 

Both  clean  and  control  are  written  in  the  ANSI/C  language,  and  compi'e  and  operate  under  ei¬ 
ther  the  System  V  (HPUX),  BSD4.2,  or  BSD4.3  versions  of  the  Unix  operating  system. 

2.2.4  Current  Status 

Version  v2.3  of  the  clean  and  control  software  and  their  corresponding  documentation  sets 
were  delivered  to  NAVOCEANO  in  June  of  1992.  On-site  software  training  was  provided  to 
NAVOCEANO  personnel  at  that  time.  The  clean/control  software  and  User’s  Guide  were  also  de¬ 
livered,  at  DMA’s  request,  to  a  NOAA  facility. 
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Items  delivered  under  this  task  include  the  software  source  code,  functional  system  and  sub¬ 
system  descriptions  |8,9j,  and  User’s  Guide  { 10}.  The  final  software  release  will  be  sent  to  DMA 
at  project  completion,  circa  1  May  93. 

2.3  Task  3:  Feasibility  Study 

In  this  task,  SRI  conducted  a  study  to  identify  data  processing  techniques  already  in  use  for 
other  scientific  applications  that  could  be  cost-effectively  adapted  to  oceanography.  Wherever 
possible,  estimates  were  made  for  the  additional  required  R&D,  implementation  time/cos.t  and 
technical  risk  for  adapting  such  techniques  to  bathymetric/hydrographic  data  processing  prob¬ 
lems.  In  the  particular  area  of  error  reduction  through  sonar  signal  processing  improvements,  a 
prototype  algorithm  was  implemented  and  tested. 

The  techniques  covered  in  the  feasibility  study  fell  roughly  into  two  categories.  High-level 
data  processing  techniques,  which  are  largely  based  on  methods  of  computer  vision  and  artificial 
intelligence,  were  considered  for  the  problems  of  bathymetric  data  fusion,  seafloor  classification/ 
provincing,  and  three-dimensional  shape  and  object  recognition.  Lower-level  processing  tech¬ 
niques,  which  are  based  on  methods  of  statistics  and  mathematical  optimization,  were  studied  for 
application  to  post  facto  survey  navigation  correction,  sonar  error  reduction,  and  interferometric 
sidescan  sonar  processing. 

A  complete  report  on  our  findings  was  given  in  “Oceanographic  Data  Processing  Feasibility 
Study”  [1],  which  was  delivered  to  DMA  in  February  1993.  The  following  sections  summarize 
our  findings  for  each  of  the  application  areas.  Recommendations  for  additional  work  in  these 
areas  are  given  in  Section  3. 

2.3.1  Ocean-Bottom  Classification/Provincing 

This  portion  of  the  study  considered  methods  for  automatically  segregating  the  ocean  floor 
into  regions  that  are  similar  in  some  regard.  Typical  provincing  criteria  include  bottom 
composition,  roughness,  slope,  or  some  particular  geological  feature.  Even  if  humans  are  still 
used  to  provide  the  final  interpretation,  an  automatic  or  semiautomatic  provincing  capability 
would  be  expected  to  save  analysis  time.  For  example,  in  geophysical  studies,  the  data  analysts’s 
job  would  be  greatly  simplified  if  regions  of  similar  topology  were  identified  in  advance  by  the 
computer.  Provincing  may  also  be  used  to  automatically  select  outlier  detection  parameters,  and 
recent  findings  by  the  DMA  Data  Integration  Working  Group  (DI  WG)  indicate  that  provincing 
will  be  required  to  delineate  seafloor  regions  and  automatically  select  appropriate  gridding  and 
interpolation  algorithms.  Provincing  capabilities  of  this  type  would  also  be  useful  to  the  ocean 
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engineering  community,  where  tasks  such  as  cable  routing  could  be  simplified  if  bathymetric- 
provinces  were  automatically  delineated  and  classified. 

The  bathymetric  classification  problem  has  similarities  to  problems  that  have  already  been 
addressed  in  other  fields.  For  example,  computer  vision  techniques  are  currently  used  to  perform 
such  tasks  as  characterizing  topographies  or  assessing  crop  yields  from  digital  terrain  maps.  In 
manufacturing  applications,  computer  vision  techniques  provide  real-time  atsembly-line  inspec¬ 
tion  and  quality  control.  Such  techniques  are  well  understood,  and  could  be  transferred  to  oceano¬ 
graphic  applications  in  a  relatively  short  time. 

2.3.2  Three-Dimensional  Shape  Recognition 

A  general  three-dimensional  shape  recognition  capability  is  a  prerequisite  to  several  automatic 
tasks,  including  macroscopic  error  detection  (e.g.,  identifying  biases  in  the  sensor  estimates),  nav¬ 
igational  hazard  detection,  and  recognition  of  topographic  features. 

As  in  the  case  of  classification  and  provincing,  a  great  deal  of  research  has  already  been  con¬ 
ducted  on  terrestrial  map  interpretation,  as  well  as  on  parts  identification  and  fault  detection  in 
manufacturing  applications.  Although  a  variety  of  algorithms  have  been  developed  to  address  the 
three-dimensional  shape  recognition  problem,  they  rely  almost  universally  on  the  differential  geo¬ 
metric  concepts  of  Gaussian  and  total  curvature.  Scale  also  plays  an  important  role  in  the  shape 
recognition  process,  indicating  that  bathymetric  shape  analysis  algorithms  will  be  required  to  ex¬ 
tract  curvature  features  over  scale. 

2.3.3  Data  Fusion 

Data  fusion  is  the  process  of  combining  information  from  multiple  sensors  to  produce  com¬ 
posite  data  that  are  either  more  accurate  or  more  reliable  than  what  can  be  derived  from  any 
sensor  alone.  In  bathymetric  applications,  the  candidate  sensor  types  include  single-  and  multi¬ 
beam  sonars,  sidescan  sonars,  and  gravity  and  magnetic  sensors.  Fusion  may  take  place  at  low  or 
high  levels  of  data  abstraction.  For  example,  depth  data  from  an  overlapping  set  of  sonar  surveys 
may  be  combined  at  the  raw  depth-estimate  level.  Techniques  for  combining  these  data  are  well- 
known  and  relatively  straightforward  to  implement;  however,  combining  data  from  sensors  that 
measure  very  different  phenomena  (e.g.,  magnetic  and  depth  data)  is  much  more  difficult.  In  most 
cases,  features  or  underlying  interpretations  must  be  extracted  from  each  sensor’s  data  stream 
prior  to  fusion. 

An  initial  practical  goal  would  be  to  develop  a  low-level  fusion  system  that  combines  two  or 
more  bathymetric  surveys  into  a  single  map  product.  Such  a  system  could  be  developed  at 
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relatively  low  cost,  and  with  little  technical  risk.  The  development  of  a  navigational  correction 
algorithm  would  be  a  prerequisite  to  this  effort.  Additional  research  will  be  required,  however, 
before  a  high-level  fusion  system  can  be  developed.  In  particular,  knowledge  representations  must 
be  identified  to  support  specific  oceanographic  fusion  problems  and  sensor  types. 

2.3.4  Error  Reduction  Through  Improved  Sonar  Signal  Processing 

The  signal  processing  and  bottom  estimation  algorithms  currently  used  in  bathymetric  sonar 
systems  are  subject  to  outlier  errors  under  a  number  of  common  conditions.  Sonar  beam  side- 
lobes,  external  interferers,  and  poor  bottom  reflectivity  can  greatly  affect  system  accuracy,  leading 
to  such  well-known  problems  as  “tunneling,”  “omega”  effects,  phantom  features,  and  isolated  out¬ 
liers.  While  some  of  these  errors  can  be  detected  in  postanalysis,  many  are  insidious  and  are  easily 
missed  by  an  analyst.  And  when  such  errors  are  detected,  depth  data  in  the  affected  region  are  un¬ 
recoverable,  and  the  area  must  be  resurveyed  or  eliminated  from  the  resulting  data  products. 

de  Moustier  and  Kleinrock  f  1 1  ]  pointed  out  that  the  signal  due  to  the  acoustic  return  from  the 
seafloor  is  generally  still  present  in  the  data,  but  a  much  stronger  signal  at  an  incorrect  apparent 
time-delay  confuses  the  depth  estimation  algorithms.  It  is  believed  that  sophisticated  depth  esti¬ 
mation  techniques  could  make  sonar  systems  much  more  robust  against  such  problems  by  com¬ 
paring  the  shape  and  spatial  spectra  of  the  seafloor  vs.  those  caused  by  various  error  effects,  as 
well  as  the  amplitude  characteristics  of  actual  and  error  signals. 

The  problem  of  robust  detection  and  estimation  has  already  been  faced  in  several  other  appli¬ 
cation  areas.  In  the  field  of  robotics,  surface  reconstruction  algorithms  have  been  developed  to  es¬ 
timate  visible  surfaces  from  range  measurements  that  are  very  similar  to  those  produced  by  bathy¬ 
metric  sonars.  Such  approaches  rely  on  a  combination  of  surface  shape  and  consistency 
constraints  to  identify  likely  surface  interpretations.  The  target  tracking  algorithms  used  in  radar 
systems  address  the  robustness  problem  through  use  of  a  multiple  hypothesis  (MH)  approach,  in 
which  several  candidate  target  tracks  (or,  in  a  bathymetric  application,  several  candidate  sea¬ 
floors)  are  developed  and  evaluated  in  parallel.  The  final  selection  is  delayed  until  after  sufficient 
supporting  data  have  been  analyzed.  The  MH  approach  is  particularly  useful  for  applications  in 
which  sensor  errors  are  not  always  immediately  obvious,  such  as  in  the  cases  of  tunneling  and 
omega  errors  in  a  muitibeam  survey. 

A  prototype  depth  estimation  procedure  was  developed  and  briefly  tested  under  this  study. 

The  prototype  combines  a  multiple  hypothesis  procedure  with  a  surface  reconstruction  algorithm 
that  was  derived  from  the  computer  vision  literature  (e.g.,  Amini  et  al.  [12,13]).  Preliminary 
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experimental  results  using  synthesized  depth  data4  have  been  encouraging,  although  additional 
tests  on  collected  sonar  data  will  be  required  to  thoroughly  prove  the  concept. 


2.3.5  Automated  Navigational  Correction 

Prior  to  the  introduction  of  the  Global  Positioning  System  (GPS),  navigational  errors  were  the 
predominant  source  of  uncertainty  in  bathymetric  surveys.  Unfortunately,  most  of  the  existing 
bathymetric  database  was  created  in  the  pre-GPS  era,  and  many  of  those  surveys  contain  ship  po¬ 
sition  errors  of  up  to  several  kilometers.  While  most  bathymetric  map  production  facilities  allow 
navigational  corrections,  such  procedures  are  manual  and  time-consuming.  In  addition,  bathy¬ 
metric  data  must  sometimes  be  discarded  when  surveys  overlap  but  do  not  precisely  align,  be¬ 
cause  there  is  no  low-level  fusion  capability  that  accounts  for  the  navigation  errors  and  the  re¬ 
sulting  disagreement  in  depths.  For  these  reasons,  the  bathymetrists  interviewed  under  this  effort 
indicated  great  interest  in  an  automatic,  postsurvey  navigational  correction  capability.  A  similar 
capability  for  verifying  (but  not  attempting  to  correct)  the  navigation  of  coastline  surveys  was  also 
requested  by  bydrographers. 

A  search  of  the  open  technical  literature  found  that  at  least  two  prior  researchers  developed 
post  facto  navigation  correction  techniques  [15-171.  In  both  techniques,  the  seafloor  topography  is 
used  as  an  alignment  feature  at  points  where  the  survey  track  passes  over  itself  or  over  an  indepen¬ 
dent  (and  possibly  more  accurate)  survey.  These  approaches  reduce  the  relative  position  errors 
that  occur  over  the  period  of  the  survey,  and  if  the  navigation  errors  are  at  least  partially  uncorre¬ 
lated,  the  absolute  position  errors  are  reduced  as  well.  More  importantly,  if  two  or  more  survey 
points  are  precisely  known  (e.g.,  through  GPS  fixes,  or  when  a  survey  track  overlaps  a  separate, 
GPS-navigated  survey),  the  adjusted  positional  accuracy  of  the  entire  track  can  approach  that  of  a 
completely  GPS-navigated  survey. 

Based  on  the  results  of  previous  research  in  this  area,  we  conclude  that  a  postsurvey  naviga¬ 
tion  system  is  practical.  However,  before  such  a  system  can  be  used  routinely  at  an  oceanographic 
map-production  facility,  several  problem  areas  will  have  to  be  addressed,  including:  robustness  to 
sensor  errors,  especially  those  resulting  in  phantom  topography  (e.g.,  omega  effects);  multiple  ap¬ 
parent  topographic  feature  alignments,  a  situation  that  is  common  in  areas  with  few  significant 
features;  and  numeric  instabilities  associated  with  solving  problems  of  this  size. 


4.  The  sonar  data  necessary  for  these  experiments  were  created  by  combining  actual  seafloor  depth  measurements 
with  a  sonar  signal  synthesizer.  The  sonar  signal  generator  was  taken  from  a  model  given  by  Okino  and  Higashi  [14], 


2.3.6  Interferometric  Sidescan  Sonar  Data  Processing 

Interferometric  sidescan  sonar  swath  bathymetry  has  the  advantages  of  low-cost,  simple  de¬ 
sign,  and  wide  swaths  that  make  it  an  attractive  alternative  to  multibeam  bathymetry  systems. 
However,  the  simplicity  of  the  ohysical  system  is  offset  by  the  complexity  of  the  processing  re¬ 
quired  to  extract  accurate  bathymetry  from  reflected  acoustic  signals.  If  interferometric  sonar  is  to 
produce  a  reliab'e  bathymetric  produ''*  that  competes  favorably  whh  multibeam  data,  larce  im¬ 
provements  in  the  processing  stream  must  occur.  These  improvements  must  be  implemented  at  all 
stages  of  the  production  process,  from  raw  data  processing  to  fusion  with  other  bathymetric  data 
sources.  In  addition  to  the  problems  shared  with  muitibeam  systems  that  we  have  addressed  else¬ 
where  in  this  report,  we  have  identified  for  future  research  several  areas  that  are  unique  to  interfer¬ 
ometric  systems,  as  described  below. 

Optimal  Phase  Recovery.  Bathymetry  is  recovered  from  phase  in  an  interferometric  system. 
Unfortunately,  phase  is  a  random  process  influenced  by  multiple  factors.  Some  attempt  should  be 
made  to  apply  optimal  estimation  techniques  to  the  phase  recovery  problem. 

Standard  Processing  Sequence.  Because  interferometric  systems  are  just  now  emerging  from 
the  research  stage,  no  production  processing  sequence  has  yet  been  developed.  A  standard  ap¬ 
proach,  drawn  from  lessons  learned  in  the  multibeam  bathymetry  field,  should  be  developed. 

Beam  Sharpening.  Sidescan  sonar  is  subject  to  the  blurring  effects  of  beam  patterns.  These  ef¬ 
fects  can  be  reduced  in  radiometric  images  by  deconvolution.  Studies  should  be  conducted  to  de¬ 
termine  if  similar  techniques  can  be  applied  to  bathymetric  data. 

SVP  Error  Reduction.  Interferometric  bathymetry  is  very  sensitive  to  changes  in  sound  ve¬ 
locity  profile  (SVP),  requiring  frequent  recalibration  on  flat  bottom.  Techniques  that  reduce  the 
frequency  of  recalibration  or  ameliorate  the  effects  of  rapidly  changing  conditions  would  greatly 
improve  the  quality  of  data. 
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3  CONCLUSIONS 


The  outlier  detection  software  was  successfully  integrated  with  the  bathymetric  production 
system  at  NAVOCEANO,  where  it  is  now  used  routinely.  The  outlier  detection  algorithm  im¬ 
provements  recommended  under  our  earlier  program  were  also  implemented,  and  a  preliminary 
evaluation  was  carried  out.  Based  on  these  tests,  the  primary  goal  of  improving  outlier  detection 
performance  near  regions  of  significant  topographic  features  appears  to  have  been  met.  An  exten¬ 
sive  performance  evaluation  using  human  experts  as  a  control,  however,  was  not  part  of  this  pro¬ 
gram. 

The  feasibility  studies  carried  out  under  this  effort  indicate  that  a  significant  body  of  ideas  and 
algorithms  is  available  from  existing  computer  vision  and  surveillance  data  analysis  systems. 
Problems  such  as  data  fusion,  data  alignment,  and  object/feature  recognition  have  already  been 
successfully  addressed  for  these  applications.  Additional  research  and  a  significant  amount  of  en¬ 
gineering  will  be  required  before  advanced  techniques  of  this  kind  can  be  incorporated  into  ocean¬ 
ographic  data  production  systems. 

To  gain  perspective  in  the  feasibility  study  on  the  interactions  and  prerequisites  for  oceano¬ 
graphic  data  processing,  a  generic  data  processing  block  diagram  was  developed  around  the  re¬ 
quirements  put  forth  by  the  DHI  program.  As  shown  in  Figure  1 1,  the  overall  system  comprises 
three  levels:  Sonar  Processing,  which  includes  depth  estimation  and  shipboard  navigation;  Data 
Conditioning  and  Assessment,  which  performs  data  editing  (outlier  removal,  etc.),  renavigation, 
and  quality  assessment;  and  Product  Generation,  which  creates  the  finished  data  products;  i.e., 
high-level  data  fusion  to  support  geophysical  sciences,  provincing/classification  for  ocean  engi¬ 
neering  and  geomorphological  studies,  shape  recognition  to  create  feature  databases  and  locate 
potential  hazards  to  navigation,  and  low-level  fusion  and  gridding  to  create  bathymetri^ydro- 
graphic  map  products.  Our  investigations  found  that  few  of  the  processing  functions  outlined  in 
Figure  1 1  currently  exist;  of  those  shown,  only  the  data  editing  and  data  gridding  blocks  have 
been  developed  and  standardized5  for  the  oceanographic  community. 

Table  3  summarizes  the  processing  technologies  that  will  be  required  to  support  the  DHI  data 
production  capabilities.  Note  that  each  of  the  desired  capabilities  will  require  the  support  of  sev¬ 
eral  processing  technologies,  and  conversely  (and  perhaps  more  importantly),  each  processing 
technology  will  also  be  useful  to  several  application  areas.  We  also  believe  that  the  development 
of  processing  technologies  will  be  expensive,  and  care  must  be  taken  to  reduce  the  redundancy  of 
future  R&D  efforts  that  involve  them. 


5.  A  standardization  effort  for  gridding  algorithms  is  currently  under  way. 
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FIGURE  1 1  GENERIC  BATHYMETRIC  DATA  REDUCTION  MODEL 
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For  these  reasons,  we  believe  that  the  most  cost-effective  approach  to  achieving  the  DHI  goals 
is  to  identify  and  develop  critical  core  processing  technologies,  and  then  to  integrate  these  into 
each  of  the  application  areas.  Each  core  technology  should  be  developed  under  a  single  R&D  pro¬ 
gram;  a  set  of  programs  could  then  address  each  of  the  specific  application  areas.  For  example,  a 
single  shape-recognition  R&D  program  could  develop  the  tools  necessary  for  macroscopic  error 
detection  and  topographic  classification.  Similarly,  a  single  navigational  correction  capability 
could  be  developed  to  support  the  several  forms  of  data  fusion,  as  well  as  to  achieve  more  precise 
surveys. 

A  final  point  regards  the  benefits  of  involving  the  end-users  in  the  system  design  process.  Es¬ 
pecially  as  DMA  begins  to  focus  on  more  sophisticated  processing  technology,  the  required  ex¬ 
pertise  in  the  diverse  areas  of  oceanography  and  processing  are  unlikely  to  be  represented  in  any 
one  individual.  To  counter  this  problem,  the  algorithm/software  experts  must  work  closely  with 
the  oceanographers  during  the  design  stages.  Furthermore,  a  prototye-and-iterate  approach  is  rec¬ 
ommended.  In  such  an  approach,  a  low-cost  software  prototype  is  delivered  to  user  facilities  prior 
to  the  final  software  product  design,  evaluated  by  experts,  and  then  design  modifications  are 
made.  Although  adding  slightly  to  the  development  cost,  this  approach  insures  that  the  software 
addresses  the  user  community’s  needs  and,  in  the  long  run,  better  meets  the  intended  goal  of  im¬ 
proved  data  analysis  efficiency. 
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